模拟退火算法

模拟退火算法 (Simulated Annealing，SA) 最早的思想是由 N. Metropolis 等人于1953年提出。1983年, S. Kirkpatrick 等成功地将退火思想引入到组合优化领域。它是基于 Monte-Carlo 迭代求解策略的一种随机寻优算法，其出发点是基于物理中固体物质的退火过程与一般组合优化问题之间的相似性。

模拟退火算法从某一较高初温出发，伴随温度参数的不断下降,结合概率突跳特性在解空间中随机寻找目标函数的全局最优解，即在局部最优解能概率性地跳出并最终趋于全局最优。

与具有陷入局部最小值的缺点的基于梯度的方法(爬山法)以及其他确定性搜索方法不同，SA 的主要优势在于其具有避免陷入局部最小值的能力。事实上，已经证明，如果足够的随机性与非常缓慢的冷却结合使用，模拟退火将收敛到其全局最优性。

模拟退火在建模比赛中最主要的应用是模型求解，特别是最优化问题（组合优化问题），如旅行商问题 ( TSP）、最大截问题(Max Cut Problem)、0-1背包问题(Zero OneKnapsack Problem)、图着色问题(Graph Colouring Problem)、调度问题(Scheduling Problem)等等。

物理中的退火

退火是指将固体加热到足够高的温度，使分子呈随机排列状态，然后逐步降温使之冷却，最后分子以低能状态排列，固体达到某种稳定状态。

加温过程：增强粒子运动，消除系统原先可能存在的非均匀态。

等温过程：对于与环境换热而温度不变的封闭系统，系统状态的自发变化总是朝自由能减少的方向进行，当自由能达到最小时，系统达到平衡。

冷却过程：使粒子热运动减弱并渐趋有序，系统能量逐渐下降，从而得到低能的晶体结构。

Boltzman 概率分布

在统计学中，玻尔兹曼分布的表达形式为 $$ F(state) \propto e^{-E/kT} $$ 上式中，$E$是能量状态，$kT$是玻尔兹曼常数和热力学温度的乘积。

系统的两种状态之间的玻尔兹曼分布比率称为玻尔兹曼因子， $$ \dfrac{F(state2) }{F(state1) } = e^{-(E_1 - E_2 )/kT} $$

从上述表达式中可以得到：

• 同一个温度，分子停留在能量小状态的概率大于停留在能量大状态的概率
• 温度越高，不同能量状态对应的概率相差越小，温度足够高时，各状态对应概率基本相同。
• 随着温度的下降，能量最低状态对应概率越来越大，温度趋于0时，其状态趋于1

劣解接受概率（Metropolis准则）

又称做状态接受函数，这里是模拟退火法这个名字的来源，模仿固体退火原理，随着温度的下降（迭代次数上升），能量逐渐稳定，即劣解的接受概率p逐渐下降，其公式为： $$ \displaystyle p=\left\{\begin{array}{ll}{1,} & {E\left(x_{new}\right) < E\left(x_{o l d}\right)} \\ {exp \left({-\dfrac{E\left(x_{n e w}\right)-E\left(x_{o l d}\right)}{T}}\right),} & {E\left(x_{n e w}\right) \geq E\left(x_{o l d}\right)}\end{array}\right. $$

从公式中可以看出这是求最小值时的分布，因为当新的解小于旧解时百分百接受。又可以看出当$E\left(x_{n e w}\right) \geq E\left(x_{o l d}\right)$时，$p$恒在[0, 1]之间。可以看到，随着时间推移，温度越低，对于劣解的接受概率越低。

在同一个温度 $T$, 选定两个能量 $E_{1}<E_{2}$, 有 $$ P\left\{\bar{E}=E_{1}\right\}-P\left\{\bar{E}=E_{2}\right\}=\frac{1}{Z(T)} \exp \left(-\frac{E_{1}}{k_{B} T}\right)\left[1-\exp \left(-\frac{E_{2}-E_{1}}{k_{B} T}\right)\right] $$ 在同一个温度, 分子停留在能量小的状态的概率比停 留在能量大的状态的概率要大。

Metropolis准则 (1953）–以概率接受新状态

固体在恒定温度下达到热平衡的过程可以用Monte Carlo方法（计算机随机模拟方法）加以模拟，虽然该方法简单, 但必须大量采样才能得到比较精确的结果, 计算量很大。

若在温度 $T$, 当前状态 $i \rightarrow$ 新状态 $j$ 若 $E_{j}<E_{i}$, 则接受 $j$ 为当前状态; 否则，若概率 $p=\exp \left[-\left(E_{j}-E_{j}\right) / k_{B} T\right]$ 大于 $[0,1)$ 区间的随机 数, 则仍接受状态 $j$ 为当前状态; 若不成立则保留状态 $i$ 为当前状态。

$$ p=\exp \left[-\left(E_{j}-E_{i}\right) / k_{B} T\right] $$

在高温下, 可接受与当前状态能量差较大的新状态; 在低温下，只接受与当前状态能量差较小的新状态。

组合优化与物理退火的相似性

$$ \begin{array}{c|c} \hline \text { 组合优化问题 } & \text { 金属物体 } \\ \hline \text { 解 } & \text { 粒子状态 } \\ \hline \text { 最优解 } & \text { 能量最低的状态 } \\ \hline \text { 设定初温 } & \text { 熔解过程 } \\ \hline \text { Metropolis抽样过程 } & \text { 等温过程 } \\ \hline \text { 控制参数的下降 } & \text { 冷却 } \\ \hline \text { 目标函数 } & \text { 能量 } \\ \hline \end{array} $$

算法描述

• 1) 随机产生一个初始解 $x_{0}$, 令 $x_{\text {best }}=x_{0}$, 并计算目标函数值 $E\left(x_{0}\right)$;
• 2) 设置初始温度 $T(0)=T 0$, 迭代次数 $i=1$;
• 3) Do while $T(i)>T_{m i n}$
  • 1) for $j=1 \sim k$
  • 2) 对当前最优解 $x_{\text {best }}$ 按照某一邻域函数, 产生一新的解 $x_{\text {new }}$ 。 计算新的目 标函数值 $E\left(x_{\text {new }}\right)$, 并计算目标函数值的增量 $\Delta E=E\left(x_{n e w}\right)-E\left(x_{\text {best }}\right)$ 。
  • 3) 如果 $\Delta E<0$, 则 $x_{\text {best }}=x_{\text {new }}$;
  • 4) 如果 $\Delta E>0$, 则 $p=\exp (-\Delta E / T(i))$; -1) 如果 $c=\operatorname{random}[0,1]<p, x_{\text {best }}=x_{\text {new }}$; 否则 $x_{\text {best }}=x_{\text {best }}$ 。
  • 5) End for
• 4) $i=i+1$;
• 5) End Do
• 6) 输出当前最优点, 计算结束。

基本步骤

优缺点

模拟退火算法的优点

• 质量高；
• 初值鲁棒性强；
• 简单、通用、易实现。 模拟退火算法的缺点
• 由于要求较高的初始温度、较慢的降温速率、较低的终止温度，以及各温度下足够多次的抽样，因此优化过程较长。

---